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ABSTRACT 

A new inverse inviscid method suitable for the design of rotating blade sections 




lying on an 


arbitrary axisymmetric stream-surface with varying streamtube width is presented. 

Given are the geometry of the axisymmetric stream-surface and the streamtube width 
variation with meridional distance, the number of blades, the inlet flow conditions, the rotational 
speed and the suction and pressure side velocity distributions as functions of the normalized arc- 
length. The flow is considered irrotational in the absolute frame of reference and compressible. The 
output of the computation is the blade section that satisfies the above data. 

The method solves the flow equations on a (<J>UP) potential function-streamfunction plane 
for the velocity modulus, W and the flow angle p; the blade section shape can then be obtained as 
part of the physical plane geometry by integrating the flow angle distribution along streamlines. The 
(<£1,111) plane is defined so that the monotonic behaviour of the potential function is guaranteed, 

even in cases with high peripheral velocities. _ 

The method is validated on a rotating turbine case and used to design new blades. To obtain 
a closed blade a set of closure conditions has been developed and refered in the paper. 
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INTRODUCTION 

The design method which is presented in this paper is developed in order to use the results 
of the meridional plane calculation and in particular the geometry of the axisymmetric flow 
streamtubes. The design method is, then, applied in order to specify the blade section shape lying 
on each axisymmetric stream-surface. The complete blade is constructed by stucking these blade 
section shapes in the span-wise sence, as desired. 

Blade design methods have already been developed in the past, for both incompressible and 
compressible flows (refs [1]-[11]). However, most of them refer to plane cascade configurations only. 
During recent years the topic of developing blade design methodologies has received particular 
attention and important contributions have been published in this framework (refs [1]*[5]). 

The aim of the present effort was to develop an inverse inviscid method supporting the blade 
optimization procedure described in reference [19] and capable to deal with the general case of an 
arbitrary rotating cascade. The method follows the work of Schmidt and Zanneti concerning 
i the equations employed. However, it formulates the problem in a different way and employs i 
l different numerical techniques as well as closure conditions, for reasons explained below. A first] 
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["version of the present method, using the classical potential function/stream-function definitions was 
I presented in reference [12], Nevertheless, problems occured, when the method was applied to high s 
speed rotating cascades because of the non-monotonic behaviour of the potential function. Recently, 
a new version has been developed, capable of overcoming this problem. This improved version 
which makes use of a more appropriate definition of the potential function/stream-function plane 
(here refered as (3>l,iy)) is presented in the present paper. 

POSITION OF THE PROBLEM AND DEVELOPMENT OF THE EQUATIONS 

A schematic representation of a peripheral cascade is given in figure 1. The aim is to com- 
pute a closed blade section, given the stream-surface geometry, the streamtube width variation with 
the meridional distance (m), the^approximate number of blades (N), the inlet stagnation conditions 
(P T1 ,T T1 ) and velocity vector (W^), the meridional position of the inlet stagnation point (mj), the 
rotational speed (oo) and the derived outlet flow angle (f^)- Assumed given, as well, are the suction 
side velocity distribution and an approximate pressure side velocity distribution versus arc length. 
The number of blades and the pressure side velocity distribution will change during the computation, 
in order to obtain a closed profile, with the constraint to alter them as little as possible. 

The flow is considered steady, inviscid, compressible subsonic and irrotational in the absolute 
frame of reference. 

The physical plane is presented in figure 2a. The equations written on an axially symmetric 
system (m,0), are: 

a) the continuity equation 


3 


* (9 R (An) W m ) + — (p R (An) WJ = 0 


m 


R^0 


( 1 ) 


b) the absolute irrotational flow equation 
1 )(RW U + cu R 2 ) ^W n 
R cm Rc)0 


= 0 


( 2 ) 


In the previous version of the method (ref.[ 12]) a transformation is performed to the (<J>,IIJ) 
plane defined as 

n x V s IP = p (An) W (3) 

V s 4> = (W + co x R) (4) 

where V s is the surface gradient operator and n the normal to the surface unit vector. 

The difference d<I> along iso-W lines is equal to 

d<I> = (W + co x R) ds = Wds + coR 2 d0 (5) 

This difference, however, is not always positive since there may exist certain high peripheral 
speed cases for which d<I> locally takes negative values and, thus, 4> is non-monotonic along 
streamlines. This fact prohibits the mapping of the physical coordinate plane to the 
potential/stream-function plane (the Jacobian of the transformation becomes zero) and thus no 
arithmetical solution is possible. To overcome this inconveniency, a new transformed plane (<I>1,1P) 
is defined in the following way 
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[" Vs 01 = W + U - (B/R) ( 7 ) | 

where B is a vector parallel to the peripheral velocity U and its modulus B is constant. The 
difference d<l>l along a streamline is, then, equal to 

d<Dl = (W + toxR-B/R) ds = Wds+(coR 2 -B)d0 (8) 

It is obvious that B can be selected is such a way that guaranties the positivity of the Jacobian of 

the transformation from the physical to the (<I>l,W)-plane. 

On the (<J>l,<I»)-plane the equations of continuity and absolute irrotational flow can be 

written in the form 

2 2 
Al(lnW)* lw + A2(lnW)* 1 + A3(lnW) 01 + A4(lnW) viny + A5(lnW) w + 

A6(lnW) w + ATClnW)^ + A8(lnW) w (lnW) w + A9 = 0 (9) 


^P 

— = Fl(W,p,R,(An)) 

^ ip 

^P 

= F2(W,p,R,(An)) 

^<I>1 


( 10 ) 


( 11 ) 


The expressions for the coefficients A1 to A9 are given in the reference [15]. In the above 
equations <I>1 and U* are the independant variables, while the velocity modulus (W) and the flow 
angle (P) are the dependent ones. Equations (10) and (11) for the flow angle are equivalent so 
during the calculation one of them may be utilized. 


THE BOUNDARY CONDITIONS ON THE (4>1,V)-PLANE 

The transformed plane-(<J>l,U<) is presented in figure 2b. The flow quantities are known at 
station (1), inlet, and the flow angle at station (2), outlet. The integral mass flux conservation 
equation, the energy conservation equation along a meridional streamline and the isentropic flow 
relations are used to calculate flow quantities at station (2). The integral mass flux equation is 
written in the form 


Rj cosPj (An) t 

e 2 w 2 = <?iW 2 

R 2 cosP2 (An ) 2 


( 12 ) 


and the energy conservation equation along with the isentropic flow relations results to the following 
expression 


p 2 WjMJj 2 W 2 -U 2 l/y-l 

-= (1 + ) 

Q 1 


(13) 


2c p T x 


2CpT! 


From these two equations the flow quantities (p 2 ,W 2 ) ma y ^ calculated for a known flow angle p 2 . 
The integral momentum equation can be written in the form 
f 2n 

r = (Wds = owds + r, = — (RiV ul -R 2 v u2 ) O 4 ) 

I N 

V • 

blade blade 
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r< 


where 


fj = Q)a) R 2 d0 


blade 


This equation relates the flow conditions at the inlet and the outlet with the circulation T which 
depends on the velocity distributions along suction and pressure side, as well as, on the blade section 
geometry. Integral Tj, depending on the blade section geometry, is not a priori known. This explains 
one of the difficulties of the inverse methodology applied to arbitrary rotating cascades. Note that 
Tj is zero only when the radius R is constant. During the computational procedure the integral T l 
is given an initial reasonable value and corrected accordingly, each time a blade section shape is 
computed. In any case the value of T must be compatible with the imposed value of the outlet flow 
angle Pj, so that, if the suction side velocity distribution (being most sensitive) must be maintained, 
the pressure side velocity distribution must be chosen to satisfy this value of T. 

Considering, again, figure 2, periodic conditions are imposed along the ((AB),(EF)) and 
((CD),(GH)) pairs of boundaries. W(<I>1) is specified along the suction and pressure side solid 
boundaries and the corresponding value of Wl is calculated from the following relation 

d<I>l = Wds + (o)R 2 -B) d9 (15) 

Consequently, differences in potential from a station v to a station p may be calculated as 


A<I>1 


V- 

v 


= Wds + 


v 


(o)R 2 -B) d9 

V 


(16) 


Moreover, the way that the (<t>l,4*)-plane was built assures that 



B 
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F 

D 
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+ A<1>1 
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(17) 


During the computational procedure, the magnitudes of A4> | ® and A<t> | j? are specified with 
the constraint to take them large enough in order to reach at AE and DH (seeiigure 2b) uniform 
conditions with sufficient accuracy. In this way, the position of the inlet and outlet of the calculation 
domain in the physical plane (positions of AE and DH in figure 2a) is not yet specified. However, 
using equation (4) along the peripheral direction one may get 


dW = p(An)WcospRd9 


(18) 


so that the corresponding stream function differences are described by the following relation at the 
inlet and the outlet stations 


AW 


v 


> 

p(An)WcospRd9 

J 


(19) 


Along the inlet and outlet stations the flow is uniform with velocities and flow angles, Wj, W 2 and 
Pj, P 2 , respectively. Consequently, if W E = W F = W G = W H = 0 is the streamfunction value characterizing 
the lower boundary, then the one characterizing the upper boundary, according to equation ( 19), is 
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A4* 


2nR, 


2TtR 0 


= p 1 W 1 cosp 1 — 
E N 


(An)j = P 2 ^ 2 c0S P 2 (An) 2 =A4* 


N 


D 

H 


( 20 ) 


The upper boundary being a streamline, *P A = 1 P B =iP c SOO- 
THE NUMERICAL INTEGRATION OF THE EQUATIONS 

Equations (9) and (10) or (11) are considered in the (<I>l,4*)-plane, which in general is non- 
orthogonal. If one considers suction and pressure side extensions of equal length in the periodic 
zones (to facilitate the application of the periodicity conditions), then the computational domain 
on the (<J>l,*P)-plane takes a trapezoidal form (see figure 2). A non uniform discretization of the 
(<I>1,4*) boundary regions was found to be efficient, permitting the stretching of the grid lines in the 
near-leading and near-trailing edge regions of the blade section, where the velocity gradients are 
large. In view of the above, the resulting grid on the (<J>l,4*)-plane, composed only of straight lines, 
is generally skewed and stretched. In order to increase the generality of the solver and the accuracy 
of the solution, avoiding at the same time complexities (such as patched grid techniques), an 
additional body-fitted coordinate transformation is performed, which maps the (<I>l,4*)-plane to an 
orthogonal (£,q)-plane with square cells (see figure 2). 

The resulting equation on W in the (£,r|)-plane is discretized by use of second-order accurate 
finite-difference/ finite volume centered schemes. The discrete equation is, then, linearized, 
transfering all non-linear terms ((lnW)^, for example) to the right hand side (fixed point algorithm). 
The resulting system of algebraic equations, which has a 9-diagonal banded, non-symmetric 
characteristic matrix, is solved iteratively using the MSIP (refs [13], [14]) method (incomplete L-U 
approximate factorization procedure). 

Once the velocity field is computed, the flow angle field is obtained integrating the ordinary 
differential equations (10) or (11), along the iso-<J>l or the iso-4* lines. A fourth order Runge-Kutta 
method is used during this step. In practice, equation (11) is first integrated along the cascade mean 
streamline and the computed p-mean streamline values are used as boundary conditions for the 
integration of equation (10) along the iso-4* 1 lines. This procedure involves only a tangential 
derivation of the flow quantities along the blade and is, thus, more accurate. If a second order 
normal derivation along the blade is used, then quadratic extrapolation procedures would be 
required, decreasing the accuracy (mainly) in the sensitive leading edge region. The above procedure 
provides the complete P(4>1,4*) field and, consequently, the blade coordinates. 

THE COMPUTATIONAL ALGORITHM 

A computational algorithm was constructed, outlined by the following steps (without 
considering conditions for section closure, which will be examined later). 

STEP 1 : The exit plane flow quantities are calculated through equations (12) and (13). A value for 
the integral T, is assumed and a velocity distribution for the pressure side compatible with the value 
of the circulation T issued from equation (14) is established. The value of constant B^is defined so 
that $1 is monotonic along streamlines. The values of the potential differences A4>1 1 “ and A4>1 1 c 
are specified. 

STEP2 : A first approximation of the (4>l,4*)-plane contour is considered and the boundary 
conditions for the velocity (through equations (16) and (19)) and the angle (utilizing plausible angle 
distributions), are specified. The interior grid points of the region (BCGF) are established using 
a simple linear procedure. In the upstream (ABFE) and downstream (CDHG) regions, the points 
on the boundaries are chosen and the grid is constructed, so that periodic conditions can be checked 
without interpolation. The complete velocity and flow angle fields are initialized making use of the 
values at the boundaries, through a linear interpolation. An initial estimate of (An) and R for each 
node is made, as well. 

STEP 3 : The coefficients Ai(i= 1,9) appearing in equation (9) are calculated. , 

[ STEP 4 : Equation (9) is solved for W(4>1,4>) using the numerical procedure and technique] 
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described in the previous section. At this point, an iterative procedure is performed involving the 
previous step, that is, updating the values of the coefficients Ai. This updating is performed, utilizing 
the values of the velocity field of the previous iteration. 

At the end of the computational procedure involved in this step, the values of W at the 
periodic boundaries will have been updated along with the complete velocity field. 

STEP 5 : The flow angle field p(<I>l,W) is computed after numerical integration of equations (10) 
and/or (11) in the manner described in the previous section. During this procedure, new angle 
values are computed at the boundaries, as well. 

STEP 6 : The blade section shape 0 = 0(m) is computed using the following geometrical relations, 
valid along a streamline 


m = 


0 = 


cospds = m(s) 
sinp 


R 


ds = 0(s) 


( 21 ) 

(22) 


Utilizing these relations, the values of m and 0 are computed along streamlines for the whole flow 
field, as well. An interpolation procedure is used in order to estimate the new set of values R(m(s)) 
and An(m(s)), which will be used, along with the updated values of the angles. 

The exit conditions are calculated at station (2), using the same procedure as in STEP1. The 
integral Tj is then computed and its new value is used to update T. The pressure side velocity 
distribution is in turn modified in order to satisfy the new value of the circulation. The B constant 
value is modified for the new geometry and velocity distribution. The boundaries and associated 
conditions can then be established for a new (<£l,ip)-plane. A new grid is thus generated on the 
(«l>l,U , )-plane, moving along W-lines and computing each time the value of <I>1 corresponding to the 
previously updated values of the velocity field. 

STEP 7 : STEPS 3 to 6 are repeated until convergence is achieved. 

As observed before, the blade section shape obtained from the above described 
computational procedure is not necessarily closed. 

RESULTS AND DISCUSSION 

To validate the method stationary and rotating cascade reconstruction test cases were 
selected. Exact cases were prefered where possible, while a direct solver was used to calculate the 
"target", velocity distribution when the later was not analytically known. Inevitably, slight inaccuracies 
in the results of the direct calculation method resulted in inaccuracies of the computed blade shape 
by the inverse method. A complete outline of the test cases utilized for the validation of the method 
are reported by Bonataki^ 15 l Results for two analytical test cases and for a radial inflow turbine are 
presented below. 

In figure 3 the Gostelowf 16 ! exac t case (incompressible flow, compressor cascade) and in 
figure 4 the Hobson^ 17 ! exact case (high Mach number, high turning angle, low pitch to chord ratio) 
are presented to demonstrate the accuracy of the method. A radial inflow turbine case^ 18 ! (strong 
variation of R(m), rotational, variation of An(m)) is presented in figure 5. In all three cases the 
presented results include the initial blade shape, the corresponding suction and pressure side velocity 
distribution and the blade shape provided by the inverse method. The typical number of grid points 
utilized for the above calculations was (78x15) and the computing time needed for the complete 
solution was 20 cpu seconds in an ALLLANT FX 80 computer. 

As a next step the method was used for the design of new profiles. Starting from an arbitrary 
suction and pressure side velocity distribution, a procedure was developed which in few iterations 
| provides a closed profile. This procedure is based on an extended investigation upon the parameters. 
Lwhich influences the blade section shapel 12 ^ 15 !, an investigation which has pointed out that the ratio' 
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[of pressure to suction side arc length and the pitch to chord ratio could control blade section j 
I closure. It was also observed that the velocity distribution near the blade section edges influences 
a small part of the blade shape near these edges, while the blade thickness is directly related to the 
mean value of the velocity distributions along both the pressure and suction sides' 121 . Using the 
information provided by the above investigation new profiles were designed. 

A rotating turbine cascade lying on a conical surface along with the "target" velocity 
distribution is presented in figure 6. 

In figure 7 a turbine blade is presented, which was used as the starting point for the design 
of a thicker blade. This new blade was obtained by increasing the level of the suction and pressure 
side velocity distributions while retaining the same inlet and outlet flow conditions. This particular 
design is quite revealing, since the "target" velocity was obtained by modifying the original one in 
such a way, so that the maximum velocity along the blade surfaces was not increased. 

The blade section shape of a radial inflow turbine with speed of rotation is presented in 
figure 8a along with the corresponding R(m) and An(m) distributions (figures 8c, 8d) and "target" 
velocity distribution (figure 8b). This is a typical case where the classical (^'P) plane definition fails 
and this is demonstrated in figures 8e,8f where the (<J>l,9J)-plane is plotted for two different values 
of the B parameter, B = 0 (the classical 4> definition) and B=l.l (the modified definition). It is 
evident that the modified definition suits better to the specific case. 

CONCLUSIONS 

A new inverse inviscid method for designing stationary or rotating, plane or axisymmetric 
cascades was presented in this paper. 

Compared with previous efforts, the new method may handle cascades rotating with high 
speed and provide closed blade shapes in few external iterations. 

The formulation and the numerics of the corresponding inverse method were discussed, in 
order to distinguish it from similar methods and reveal its relative merits. 

Finally, some calculation results were presented, to certify the accuracy and the capabilities 
of the present effort. 
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Fig.l A schematic representation 
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Fig. 2 The schematic representation of a 
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Fig. 3 The Goste.low cascade 
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Fig. 5 A rotating turbine case 
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Fig. 6 Design of a new rotating cascade 
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